####################################################################################

#What You See and What You Get
#This file conducts the analysis and produces the Tables for Appendix L on MCMV allocation

####################################################################################

rm(list=ls())
library(tidyverse)
library(rdrobust)
library(modelsummary)
library(rddensity)
library(estimatr)
library(xtable)
options(scipen=999)

source("Code/functions.R")


load("Data/empresas_g_noimppt.Rda")
load("Data/empresas_gpp.Rda")
load("Data/lai2.Rda")

######################################################
############### RDD  
######################################################
#Extra help functions from here: https://www.appsloveworld.com/r/100/587/print-tables-with-3-regression-output-models-from-rdrobust

tidy.rdrobust <- function(model, ...) {
  ret <- data.frame(
    term = row.names(model$coef),
    estimate = model$coef[, 1],
    std.error = model$se[, 1],
    p.value = model$pv[, 1]
  )
  row.names(ret) <- NULL
  ret
}

glance.rdrobust <- function(model, ...) {
  ret <- data.frame(
    Kernel = model$kernel,
    Bandwidth = model$bwselect
  )
  ret
}


# RDD PT

#Some recoding for RDD data
empresas_g_noimp <- ungroup(empresas_g_noimp)
empresas_gpp <- ungroup(empresas_gpp)

#Recoding eletion year
ANO_ELEICAO2 <- factor(empresas_g_noimp$ANO_ELEICAO)
ANO_ELEICAO.d = model.matrix(~ANO_ELEICAO2+0)
ANO_ELEICAO.d = ANO_ELEICAO.d[,1]

tpt <- rdrobust(y = empresas_g_noimp$valor_contratadodpcl2,
                x = empresas_g_noimp$vote_margin_share, 
                        covs = ANO_ELEICAO.d, cluster = empresas_g_noimp$TSEcod)

#RDD PP

#Recoding eletion year
ANO_ELEICAO2 <- factor(empresas_gpp$ANO_ELEICAO)
ANO_ELEICAO.d = model.matrix(~ANO_ELEICAO2+0)
ANO_ELEICAO.d = ANO_ELEICAO.d[,1]

tpp <- rdrobust(y = empresas_gpp$valor_contratadodpcl2, x = empresas_gpp$vote_margin_share, 
                        covs = ANO_ELEICAO.d,
                        cluster = empresas_gpp$TSEcod)

models <- list(tpt, tpp)

modelsummary(models,
             stars = TRUE, 
             output = "Tables/tab-l4.tex")


#Balance recoding
empresas_g_noimp <- empresas_g_noimp %>% mutate(p_13.1 = p_13/p_val,
                                                p_45.1 =  p_45/p_val,
                                                f.nom_13.1 = f.nom_13/f.nom_total,
                                                f.nom_45.1 = f.nom_45/f.nom_total,
                                                g_13.1 = g_13/g_val,
                                                g_45.1 = g_45/g_val, 
                                                e.nom_13.1 = e.nom_13/e.nom_total,
                                                e.nom_45.1 = e.nom_45/e.nom_total,
                                                Comparecimento1t.1 = Comparecimento1t/pop_v2)


outcomes1 <- list(
  empresas_g_noimp$X2000_incomepercapita,
  empresas_g_noimp$X2000_doctor,
  empresas_g_noimp$X2000_idh_education,
  empresas_g_noimp$X2000_idh_income,
  empresas_g_noimp$X2000_idh_longevity,
  empresas_g_noimp$X2000_illiterate,
  empresas_g_noimp$X2000_infant,
  empresas_g_noimp$X2000_pop, 
  empresas_g_noimp$X2000_poverty,
  empresas_g_noimp$p_13.1, 
  empresas_g_noimp$p_45.1,
  empresas_g_noimp$f.nom_13.1, 
  empresas_g_noimp$f.nom_45.1,
  empresas_g_noimp$g_13.1, 
  empresas_g_noimp$g_45.1,
  empresas_g_noimp$e.nom_13.1,
  empresas_g_noimp$e.nom_45.1,
  empresas_g_noimp$Comparecimento1t.1)

#Minister's Alignment
empresas_gpp <- empresas_gpp %>% mutate(p_13.1 = p_13/p_val,
                                        p_45.1 =  p_45/p_val,
                                        f.nom_13.1 = f.nom_13/f.nom_total,
                                        f.nom_45.1 = f.nom_45/f.nom_total,
                                        g_13.1 = g_13/g_val,
                                        g_45.1 = g_45/g_val, 
                                        e.nom_13.1 = e.nom_13/e.nom_total,
                                        e.nom_45.1 = e.nom_45/e.nom_total,
                                        Comparecimento1t.1 = Comparecimento1t/pop_v2)


outcomes2 <- list(
  empresas_gpp$X2000_incomepercapita,
  empresas_gpp$X2000_doctor,
  empresas_gpp$X2000_idh_education,
  empresas_gpp$X2000_idh_income,
  empresas_gpp$X2000_idh_longevity,
  empresas_gpp$X2000_illiterate,
  empresas_gpp$X2000_infant,
  empresas_gpp$X2000_pop, 
  empresas_gpp$X2000_poverty,
  empresas_gpp$p_13.1, 
  empresas_gpp$p_45.1,
  empresas_gpp$f.nom_13.1, 
  empresas_gpp$f.nom_45.1,
  empresas_gpp$g_13.1, 
  empresas_gpp$g_45.1,
  empresas_gpp$e.nom_13.1,
  empresas_gpp$e.nom_45.1,
  empresas_gpp$Comparecimento1t.1)

names <- c("Income per capita", 
           "Doctors per thousand pop.",
           "Education (IDH)",
           "Income (IDH)", 
           "Longevity (IDH)", 
           "Illiteracy rate", 
           "Infant mortality", 
           "Population",
           "Poverty rate", 
           "Vote for Lula",
           "Vote for FHC", 
           "Vote for PT (fed. dep.)",
           "Votes for PSDB (fed. dep.)",
           "Vote for PT (governor)",
           "Vote for PSDB (governor)",
           "Votes for PT (state dep.)",
           "Votes for PSBD (state dep.)",
           "Turnout")

#Balance President and Minister's Alignment
poly_balance_table <- rdd.table.balance.poly(runvar = empresas_g_noimp$vote_margin_share, outcome = outcomes1) 
poly_balance_table <- do.call(rbind, poly_balance_table)
rownames(poly_balance_table) <- names

print.xtable(xtable(poly_balance_table), type = "latex", 
             file = "tables/tb-l1.tex")

poly_balance_table <- rdd.table.balance.poly(runvar = empresas_gpp$vote_margin_share, outcome = outcomes2) 
poly_balance_table <- do.call(rbind, poly_balance_table)
rownames(poly_balance_table) <- names

print.xtable(xtable(poly_balance_table), type = "latex", 
             file = "tables/tb-l2.tex")

#Manipulation tests

#President

rdd <- rddensity(X = empresas_g_noimp$vote_margin_share)
plot1 <- rdplotdensity(rdd, empresas_g_noimp$vote_margin_share)
summary(rdd)
plot1

empresas_g_08 <- empresas_g_noimp %>% filter(ANO_ELEICAO == 2008)
empresas_g_12 <- empresas_g_noimp %>% filter(ANO_ELEICAO == 2012)

rdd_08 <- rddensity(X = empresas_g_08$vote_margin_share)
rdd_12 <- rddensity(X = empresas_g_12$vote_margin_share)
summary(rdd_08)
summary(rdd_12)

pdf("figures/fig-l1a.pdf",
    width=20, height=20)
rdplotdensity(rdd_08, empresas_g_08$vote_margin_share)
dev.off()

pdf("figures/fig-l1b.pdf",
    width=20, height=20)
rdplotdensity(rdd_12, empresas_g_12$vote_margin_share)
dev.off()

#Minister

rdd <- rddensity(X = empresas_gpp$vote_margin_share)
plot1 <- rdplotdensity(rdd, empresas_gpp$vote_margin_share)
summary(rdd)
plot1

empresas_g_08 <- empresas_gpp %>% filter(ANO_ELEICAO == 2008)
empresas_g_12 <- empresas_gpp %>% filter(ANO_ELEICAO == 2012)

rdd_08 <- rddensity(X = empresas_g_08$vote_margin_share)
rdd_12 <- rddensity(X = empresas_g_12$vote_margin_share)
summary(rdd_08)
summary(rdd_12)


pdf("figures/fig-l1c.pdf",
    width=20, height=20)
rdplotdensity(rdd_12, empresas_g_12$vote_margin_share)
dev.off()

pdf("figures/fig-l1d.pdf",
    width=20, height=20)
rdplotdensity(rdd_08, empresas_g_08$vote_margin_share)
dev.off()

#SII.1 Timing of Funding 

m1 <- lm_robust(mean_dist ~ pt, fixed_effects = ~cycle + TSEcod, 
                clusters = TSEcod, se_type = "stata", data = lai2)

m2 <- lm_robust(mean_dist ~ pp, fixed_effects = ~cycle + TSEcod, 
                clusters = TSEcod, se_type = "stata", data = lai2)

m3 <- lm_robust(mean_dist ~ psdb, fixed_effects = ~cycle + TSEcod, 
                clusters = TSEcod, se_type = "stata", data = lai2)

print.xtable(xtable(rbind(tidy(m1), tidy(m2), tidy(m3))), type = "latex", 
             file = "tables/tb-l3.tex")
